COMMENT ⊗   VALID 00007 PAGES 
RECORD PAGE   DESCRIPTION
 00001 00001
 00002 00002	SUBR(SQRT)X ------------------------------------------------------
 00004 00003	INTERN SIN,COS---------------------------------------------------
 00006 00004	SUBR(ATAN)--------------------------------------------------------
 00010 00005	SUBR(ASIN)--------------------------------------------------------
 00012 00006	SUBR(LOG)---------------------------------------------------------
 00013 00007	SUBR(REALIN)
 00016 ENDMK
⊗;
SUBR(SQRT)X ------------------------------------------------------
BEGIN SQRT;MODIFIED OLDE LIB40 SQUARE ROOT.
	A←←0 ↔ B←←1 ↔ C←←2
	LACM B,ARG1↔JUMPE B,POP1J.↔PUSH P,2

;LET X=F*(2↑2B) WHERE 0.25<F<1.00 THEN SQRT(X)=SQRT(F)*(2↑B).
	ASHC B,-=27↔SUBI B,201	;GET EXPONENT IN B, FRACTION IN C.
	ROT B,-1		;CUT EXP IN HALF, SAVE ODD BIT
	DAP B,L↔LSH B,-=35	;USE THAT ODD BIT.
	ASH C,-10↔FSC C,177(B)	;0.25 < FRACTION < 1.00

;LINEAR APPROXIMATION TO SQRT(F).
	DAC C,A
	FMP C,[0.8125↔0.578125](B)
	FAD C,[0.302734↔0.421875](B)

;TWO ITERATIONS OF NEWTON'S METHOD.
	LAC B,A
	FDV B,C↔FAD C,B↔FSC C,-1
	FDV A,C↔FADR A,C
     L: FSC A,0↔LAC 1,A↔POP P,2
	POP1J↔LIT
BEND SQRT; BGB 28 DECEMBER 1972 ----------------------------------
INTERN SIN,COS;---------------------------------------------------
BEGIN SINCOS;MODIFIED OLDE LIB40 SINE & COSINE - BGB.
	A←←1 ↔ B←←2 ↔ C←←3
↑COS:	SKIPA A,ARG1
↑SIN:	SKIPA A,ARG1
	FADR  A,HALFPI			;COS(X) = SIN(X+π/2).
	MOVM B,A↔CAMG B,[17B5]↔POP1J	;FOR SMALL X, SIN(X)=X.

;B ← (ABS(X)MODULO 2π)/HALFPI
;C ← QUADRANT 0, 1, 2 OR 3.
	FDVR B,HALFPI
	LAC C,B↔FIX C,233000
	CAILE C,3↔GO[
	TRZ C,3↔FSC C,233
	FSBR B,C↔GO .-3]		;MODULO 2π.
	GO .+1(C)↔GO .+4↔JFCL↔GO[
	FSBRI B,(2.0)↔MOVNS B↔GO .+2]	;SIN(X+π)=SIN(-X)
	FSBRI B,(4.0)			;SIN(X+2π)=SIN(X)
	SKIPGE A↔MOVNS	B		;SIN(-X) = -SIN(X).

;FOR -1 ≤ B ≤ +1 REPRESENTING -π/2 ≤ X ≤ +π/2,
;COMPUTE SINE(X) APPROXIMATION BY TAYLOR SERIES.
	DAC B,C↔FMPR B,B	
	LAC A,[164475536722]↔FMP A,B
	FAD A,[606315546346]↔FMP A,B
	FAD A,[175506321276]↔FMP A,B
	FAD A,[577265210372]↔FMP A,B
	FAD A,HALFPI↔FMPR A,C↔POP1J
	LIT
BEND;-------------------------------------------------------------
INTERN HALFPI,PI,TWOPI
	HALFPI:	201622077325 ;PI/2
	PI:	202622077325 ;PI
	TWOPI:	203622077325 ;2*PI
SUBR(ATAN)--------------------------------------------------------
BEGIN ATAN
;ATAN(X) = X*(B0+A1 / (Z+B1-A2 / (Z+B2-A3 / (Z+B3))) )
;WHERE Z=X↑2, IF 0<X<=1
;IF X>1, THEN ATAN(X) = PI/2 - ATAN(1/X)
;IF X>1, THEN RH(D) =-1, AND LH(D) = -SGN(X)
;IF X<1, THEN RH(D) = 0, AND LH(D) =  SGN(X)
	A←←1 ↔ B←←2 ↔ C←←3 ↔ D←←4 ↔ E←←5
	LAC	A,ARG1		;PICK UP THE ARGUMENT IN A
ATAN1:	LACM	B, A		;GET ABSF OF ARGUMENT
	CAMG	B, A1		;IF X<2↑-33, THEN RETURN WITH...
	POP1J		;ATAN(X) = X
	HLLO	D, A		;SAVE SIGN, SET RH(D) = -1
	CAML	B, A2		;IF A>2↑33, THEN RETURN WITH
	GO[LAC A,HALFPI ↔POP1J];	ATAN(X) = PI/2
	MOVSI	C, 201400	;FORM 1.0 IN C
	CAMG	B, C		;IS ABSF(X)>1.0?
	TRZA	D, -1		;IF B ≤ 1.0, THEN RH(D) = 0
	FDVM	C, B		;B IS REPLACED BY 1.0/B
	TLC	D, (D)		;XOR SIGN WITH > 1.0 INDICATOR

	DAC B,E↔FMP B,B
	LAC C,B↔FAD C,KB3↔LAC A,KA3↔FDVM A,C
	FAD C,B↔FAD C,KB2↔LAC A,KA2↔FDVM A,C
	FAD C,B↔FAD C,KB1↔LAC A,KA1↔FDV  A,C
	FAD A,KB0↔FMP A,E

	TRNE	D, -1		;CHECK > 1.0 INDICATOR
	FSB	A, HALFPI		;ATAN(A) = -(ATAN(1/A)-PI/2)
	SKIPGE	D		;LH(D) = -SGN(B) IF B>1.0
	MOVNS A		;NEGATE ANSWER
	POP1J		;EXIT
A1:	145000000000		;2↑-33
A2:	233000000000		;2↑33

KB0:	176545543401		;0.1746554388
KB1:	203660615617		;6.762139240
KB2:	202650373270		;3.316335425
KB3:	201562663021		;1.448631538

KA1:	202732621643		;3.709256262
KA2:	574071125540		;-7.106760045
KA3:	600360700773		;-0.2647686202
BEND ATAN;--------------------------------------------------------
SUBR(ATAN2)-------------------------------------------------------
BEGIN	ATAN2

; OMEGA ← ATAN2(Y,X).
	Y←←1 ↔ X←←2
	LACM Y,ARG2↔LACM X,ARG1
	CAML Y,X↔GO L1

;HORIZONTAL TO π/2; ABS(Y) < ABS(X).
	LAC  Y,ARG2↔FDVR Y,ARG1
	PUSH 17,Y↔PUSHJ 17,ATAN		;ARCTAN(Y/X)
	SKIPL ARG1↔POP2J		;1ST & 2ND QUADRANTS.
	JUMPGE Y,[
	FSBR Y,PI↔POP2J]		;3RD QUADRANT.
	FADR Y,PI↔POP2J			;2ND QUADRANT.

;VERTICAL TO π/2; ABS(X) < ABS(Y).
L1:	LACN X,ARG1↔FDVR X,ARG2
	PUSH 17,X↔PUSHJ 17,ATAN		;ARCTAN(X/Y)
	SKIPG ARG2↔GO[
	FSB Y,HALFPI↔POP2J]
	FADR Y,HALFPI
	POP2J

BEND ATAN2;-------------------------------------------------------
SUBR(ASIN)--------------------------------------------------------
BEGIN ASIN
;ASIN(X)=ATAN(X/SQRT(1-X↑2)).
;GIVEN -1 ≤ X ≤ +1 RETURN -π/2 ≤ ASIN(X) ≤ +π/2.
	A←1 ↔ B←2
	LACN A,ARG1↔FMPR A,ARG1↔FADRI A,(1.0)
	JUMPE A,[		;WAS X EITHER -1.0 OR 1.0?
		LAC A,HALFPI
		SKIPGE ARG1
		MOVNS A↔POP1J]
	PUSH 17,A↔PUSHJ 17,SQRT
	LAC B,ARG1↔FDVR B,1↔DAC B,ARG1	;CALCULATE X/SQRT(1-X↑2)
	GO ATAN			;CALCULATE ATAN(SQRT(1-X↑2))
BEND;-------------------------------------------------------------
SUBR(ACOS)--------------------------------------------------------
;ACOS(X)= π/2 - ASIN(X).
;GIVEN -1 ≤ X ≤ +1 RETURN 0 ≤ ACOS(X) ≤ +π.
	PUSH 17,ARG1↔PUSHJ 17,ASIN
	MOVNS 1↔FADR 1,HALFPI↔POP1J
;-----------------------------------------------------------------
SUBR(LOG)---------------------------------------------------------
BEGIN LOG
	MOVM ARG1↔SKIPE 1,0↔CAMN 0,[1.0]↔POP1J
	ASHC 0,-33↔ADDI 0,211000↔MOVSM 0,TMP1#
	MOVSI 0,(-128.5)↔FADM 0,TMP1
	ASH 1,-10↔TLC 1,200000↔FAD 1,[-0.70710678]
	LAC 0,1↔FAD 0,[1.4142135]↔FDV 1,0
	DAC 1,TMP2#↔FMP 1,1
	LAC 0,[0.59897864]↔FMP 0,1
	FAD 0,[0.96147063]↔FMP 0,1
	FAD 0,[2.88539120]↔FMP 0,TMP2↔FAD 0,TMP1
	FMP 0,[0.69314718]↔LAC 1,0↔POP1J
	LIT↔VAR
BEND;-------------------------------------------------------------
SUBR(REALIN)
BEGIN REALIN;
;<EXPR>		::= <EXPR>+<TERM>|<EXPR>-<TERM>|<TERM>
;<TERM>		::= <TERM>*<PRIMARY>|<TERM>/<PRIMARY>|<PRIMARY>
;<PRIMARY>	::= -<PRIMARY>|(<EXPR>)|π|<REAL NUMBER>
	CALL(TERM)
	CAIN 1,"+"
	GO [ PUSH P,0
	     CALL(TERM)
	     FADR 0,(P)
  	     SUB P,[XWD 1,1]
	     GO REALIN+1 ]
	CAIN 1,"-"
	GO [ PUSH P,0
	     CALL(TERM)
	     MOVN 0,0
	     FADR 0,(P)
  	     SUB P,[XWD 1,1]
	     GO REALIN+1 ]
	CAIN 1,15
	CALL(GETCHL)
	POP0J
	POP0J
TERM:	CALL(PRIMARY)
TERM2:	CAIN 1,"*"
	GO [ PUSH P,0
	     CALL(PRIMARY)
	     FMPR 0,(P)
  	     SUB P,[XWD 1,1]
	     GO TERM2 ]
	CAIN 1,"/"
	GO [ PUSH P,0
	     CALL(PRIMARY)
	     EXCH 0,(P)
	     FDVR 0,(P)
  	     SUB P,[XWD 1,1]
	     GO TERM2 ]
	POP0J
;BEGIN REALIN	; INPUT SMALL REAL NUMBER - BGB - 16 DEC 1972
;AC-0 INTEGER ACCUMULATION.	AC-0 RETURNS REAL NUMBER.
;AC-1 CHARACTER.		AC-1 RETURNS BREAK CHARACTER.
;AC-2 COUNTER OF DIGITS TO RIGHT OF DECIMAL POINT PLUS ONE.
;AC-3 MINUS SIGN FLAG.
	EXTERNAL GETCHL
PRIMARY:SETZ↔SETZB 2,3
L0:	CALL(GETCHL)
	CAIN 1," "↔GO .-2
	CAIN 1,"-"↔GO[SETCMM 3↔GO L0]
	CAIN 1,"π"↔GO[MOVE 0,PI
	      GETRET: CALL(GETCHL)↔GO L3]
	CAIN 1,"("↔GO[PUSH P,3↔CALL(REALIN)↔POP P,3
		      CAIN 1,")"↔GO GETRET
		      OUTSTR[ASCIZ/WARNING: MISSING ')'
/]↔		      POP0J]
	SKIPA
L1:	CALL(GETCHL)
	CAIE 1,"."↔GO .+3↔JUMPN 2,L2↔AOJA 2,L1
	CAIL 1,"0"↔CAILE 1,"9"↔GO L2
	JUMPN 2,[CAILE 2,4↔GO L1↔AOJA 2,.+1]
	ANDI 1,17↔IMULI =10↔ADD 1↔GO L1
L2:	FLOAT↔SOSLE 2↔FDVR[1.0↔10.0↔100.0↔1000.0↔10000.0](2)
L3:	SKIPE 3↔MOVNS↔POP0J
BEND REALIN;12/16/72(BGB),14-MAR-73(TVR)-----------------------------
END;